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1. Introduction 

Future generations of ultra high bypass-ratio jet engines will require far higher pressure 
ratios and operating temperatures than those of current engines. For the foreseeable future, 
engine materials will not be able to withstand the high temperatures without some form of 
cooling. In particular the turbine blades, which are under high thermal as well as mechanical 
loads, must be cooled (Taylor, 1980, Suo, 1978 and Snyder and Roelke, 1990). Cooling of 
turbine blades is achieved by bleeding air from the compressor stage of the engine through 
complicated internal passages in the turbine blades (internal cooling, including jet- 
impingement cooling) and by bleeding small amounts of air into the boundary layer of the 
external flow through small discrete holes on the surface of the blade (film cooling and 
transpiration cooling). The cooling must be done using a minimum amount of air or any 
increases in efficiency gained through higher operating temperature will be lost due to added 
load on the compressor stage. 

Turbine cooling schemes have traditionally been based on extensive empirical data bases, 
quasi-one-dimensional computational fluid dynamics (CFD) analysis, and trial and error. With 
improved capabilities of CFD, these traditional methods can be augmented by full three- 
dimensional simulations of the coolant flow to predict in detail the heat transfer and metal 
temperatures. Several aspects of turbine coolant flows make such application of CFD difficult, 
thus a highly effective CFD methodology must be used. First, high resolution of the flow field 
is required to attain the needed accuracy for heat transfer predictions, making highly efficient 
flow solvers essential for such computations. Second, the geometries of the flow passages are 
complicated but must be modeled accurately in order to capture all important details of the 
flow. This makes grid generation and grid quality important issues. Finally, since coolant 
flows are turbulent and separated the effects of turbulence must be modeled with a low 
Reynolds number turbulence model to accurately predict details of heat transfer. 

The overall objective of our ongoing research is to develop a CFD methodology that can 
be used effectively to design and evaluate turbine cooling schemes. Grid generation for 
complicated geometries such as coolant passages, is currently an active area of research. In 
general, grid systems for complicated geometries are classified as block-structured, 
unstructured or hybrid. Of those, unstructured. grids offer the greatest flexibility for modeling 
of complex geometries and the generation of unstructured grids is largely automatic. In 
contrast, fully continuous block-structured grids, where all grid lines are at least c ! continuous 
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across block faces (here referred to as multi-block grids), are more difficult to generate but are 
the most suitable for simulations of viscous flows. In addition, flow solvers for structured 
grids typically require less memory than those for unstructured grids, and can take full 
advantage of various convergence acceleration schemes (e.g. multi grid) and fast solvers for 
implicit discretization (e.g., line Gauss-Seidel, approximate LU and ADI schemes). In this 
study, we use semi-automatically generated multiblock grids (i.e., shape of blocks is 
automatically determined but grid-topology or block-structure needs to be specified 
beforehand). 

Turbulence models used in simulations of internal flows in complicated geometries must 
be able to model flows involving separation and adverse pressure gradients. One such model is 
the k-co model of Wilcox( 1994a and 1994b). This model has several desirable features. One 
important feature is that it does not require distance to a nearest wall as a parameter. Second, 
the low Reynolds number version of the model can be used to model transition (Wilcox, 
1994b). Finally, as both k and w are well behaved numerically, stiffness associated with low- 
Reynolds number k-e turbulence models is eliminated. In addition work by Chima (1996) 
shows the model to be useful for predicting heat transfer over turbine blades. 

The present code, Glenn-HT, has been used for many types of flowfields including turbine 
blades, tip clearance flows, film cooling holes, and coolant channels. For the present 
discussion, flow and heat transfer results for the following geometries will be addressed: (1) 
rectangular ducts with a 180 degree turn, Rigby et al. (1996a), (2) a straight channel with ribs 
and bleed holes, Rigby et al. (1997), and (3) a rotating channel with ribs and a 180 degree 
turn, Rigby (1998). 

There is extensive literature available on the study of rib roughened surfaces. Some 
examples of experimental work include Boyle (1984), Chen et al. (1996), Hibbs et al. (1997), 
Rau et al. (1996), Shen et al. (1994), and Taslim et al. (1995). Some examples of numerical 
work in this area includes Stephens et al. (1995a, 1995b), Stephens et al. (1996), and Stephens 
and Shih (1997). 

For rib roughened channels with rotation, the work of Park (1996), and Park et al. (1997) 
offer local mass transfer distributions. Regional measurements of heat transfer are available in 
Wagner et al. (1992). 
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2. Overview of Glenn-HT Code (excerpts from Glenn-HT Manual) 

2.1. Introduction 

This code is a general purpose flow solver, designed for simulations of flows in 
complicated geometries. The code is based on the TRAF3D code, an efficient computer code 
designed for simulations of flows in turbine cascades (Amone et al. 1991). The Glenn-HT 
code employs the full compressible Navier-Stokes equations. It uses a multi-stage Runge- 
Kutta scheme to march in pseudo time. The code utilizes multi-grid and implicit residual 
smoothing to accelerate convergence to steady state. Convective and diffusive fluxes are 
computed using central differencing. Artificial dissipation is added to prevent odd-even 
decoupling. The discretization is formally second order accurate. To handle complex 
geometries, the code uses contiguous multiblock grid systems but has the added capability of 
handling grids with non-contiguous grid lines across branch cuts. For contiguous systems, all 
internal boundaries are conservative. The Glenn-HT code was described in detail by 
Steinthorsson et al. (1993). Note that the had been referred to as the TRAF3D.MB code until 
the name was recently changed. Some aspects of the formulation used in the code are the same 
as those described by Amone et al. (1991). For the present computations the code was fitted 
with the low Reynolds number k - to model of Wilcox (1994b). 

2.2. Governing Equations 

2.2.1. Navier-Stokes Equations 

The compressible Navier-Stokes equations for a Newtonian fluid can be written as follows 
(see e.g.. Panton, 1984): 

Continuity equation. 


^ + V.(p») = 0, 

momentum equation, 


(2.1) 


•^-(pu) + V*(pu ® u) = -Vp + V* X, (2.2) 

ot 

and total energy equation, 

J^(pe) + V»(pn(e + p/ p)) = V»£ + V«(Tn), (2.3) 
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where p is the density, m=(u,v,w) t is the velocity, p is the thermodynamic (static) pressure, 
T is the deviatoric viscous stress tensor, e is the total energy per unit mass, and q is heat flux 
due to diffusion. The last three of these quantities are given as 

t = - ^p(V«u)7 + ^p(V ® u + (V 0 u) T ), (2.4) 

where / is the 3x3 identity matrix. 


e = u + -( u • u ) , 


(2.5) 


and 


q - - kVT , (2.6) 

where p is the molecular (dynamic) viscosity, k is the thermal diffusivity, u is the internal 
energy of the gas and T is the (absolute) temperature. 

In the above equations, the operator ® is defined in terms of matrix multiples as 

u <8> v = («(v r )) r (2.7) 

where u and v are three-component column vectors. Correspondingly, the operator V ® is 
defined as 

V®v = (V(v r )) r (2.8) 

where v is a three-component column vector. The divergence operator, V* , is defined as 

V»v = (V r )v (2.9a) 

when v is a column vector and as 

V«v = ((V T )v) F (2.9b) 

when v is a three-column matrix. 

The gas is treated as ideal gas, i.e., du = c v dT , where c v is the constant- volume specific 
heat. In addition, the gas is assumed to be calorically perfect so that c v is constant, which is 
reasonably well valid for flows in turbomachinery-excluding, perhaps, flows in the 
combustor. Thus the total energy can be written as 
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( 2 . 10 ) 


e = c v T + |(m •«) 

where the reference internal energy has been taken to be zero at a temperature of absolute 
zero. 

For ideal gas, the equation of state can be written as 


and the speed of sound is 


P = P ‘KT 


( 2 . 11 ) 


c = Jyp / p (2.12) 

where is the gas constant for the particular gas mixture in use and y is the ratio of 
specific heats. 

The dynamic viscosity, p, of the gas is assumed to vary with temperature as 

£ - 

where p 0 is reference viscosity at temperature T 0 . Thermal conductivity of the gas is then 
taken as 



( 2 . 14 ) 


where c p is the constant-pressure specific heat and Pr is the Prandtl number. 

2.2.2. Formulation for Rotating Reference Frame 
To facilitate simulation of flows in rotating machinery, the governing equations are 
transformed to a coordinate system in a rotating reference frames. The advantage of this 
approach is that the grid system for the rotating geometry is stationary in the moving 
coordinate system, which leads to savings in computational effort and makes it possible to 
compute steady-state flows in the rotating frame of reference. The disadvantage is that source 
terms arise in the governing equations. 

Several different formulations have been proposed for simulating flows in rotating 
reference frames. Here, the formulation proposed by Chima and Yokota (1988) is used. In this 
formulation the governing equations are mapped to coordinate system moving with the 
rotating reference frame. The velocity vectors used as dependent variables are, however, the 
absolute velocity. Unlike the formulation of Chima and Yokota, the total energy is the absolute 
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total energy, i.e., based on absolute velocities, rather than rothalpy. Additional details can be 
found in Chima and Yokota (1988) and in Steinthorsson (1991) and Steinthorsson et al. 
(1991). 

2.3. Turbulence Model 

When working with complicated geometries it is advantageous to use a set of equations 
describing the turbulence that does not require the computation of the dimensionless distance 
to the wall y + . The boundaries between adjacent blocks should be free to cut across boundary 
layers and regions of high shear. Having to carry information on solid walls and dealing with 
comers requires communication of much information that is quite cumbersome and time 
consuming both in terms of programming and CPU time. 

The k-ca turbulence model developed by Wilcox (1994a, 1994b) satisfies our requirements. 
Modifications by Menter (1992,1993) improved the robustness of the model. Chima (1996) 
incorporated some of the latter modifications to the turbulence model and presented some 
applications of this model in the context of a Navier-Stokes solver. In fact it is the three- 
dimensional variation to the formulation adapted by Chima that has been utilized in this paper. 
Chima has shown the model to possess very good convergence properties. He also showed 
that the model performs well in predicting the rate of heat transfer from a simulated flat plate 
and turbine blades under various conditions. Below we present the dimensionless equations 
describing the turbulence in tensor notation. 

(P*,)_, + (PW +4//),j = (15) 

q.j = j=l>3 (16) 

where Sj=k and s 2 =(0 also |x,=a* . 


The source terms, P,-, of equation (1) are defined as 


p = 


Re' 1 -U,4S2 2 -^2(V. V) 


pa[\x*4£2 2 - |co( V • V)J j 
where Q is the magnitude of the vorticity. The destruction terms, are given by 


(17) 


D = 


P*pu>* 

.P P D) 2 J 


(18) 
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The coefficients appearing in the model are 

o=2.0, (3=3/40, (3*=0.09Fp, a=(5/9)(F a /F |I ), and a *=F [X , where 




5. *fr 4 
18 Ro 


Re 7 


Re n 


a n + ■ 


a A + - 


1 + 


Re T 4 


^a = ' 


1 +■ 


7?e T 




1 +■ 




( 19 ) 


Above ao=0.1, oto*=0.025,Fp=8, R^O.21, F k =6 and Re T = p*/nco. 

2.4. Boundary Conditions 

The types of boundary conditions encountered in solving the problem at hand are as 
follows: 

1) Inlet: The inlet boundary condition for subsonic flows is treated by specifying the 
relative total temperature and relative total pressure as well as the inlet angle profiles. The 
outgoing Riemann invariant is extrapolated to the inlet from within. The total temperature and 
total pressure profiles are chosen to produce specified velocity and temperature profiles. 

2) Exit: At the exit boundary, for subsonic flow, the pressure is specified while all other 
conditions are extrapolated from within. 

3) Walls: At walls, the normal pressure gradient is set to zero, either the temperature or 
heat flux is specified, and the no-slip condition is enforced. The density and total energy are 
computed from the pressure and the temperature. The boundary conditions for the turbulence 
quantities are k = 0 and 


0 ) = Sr— 

% 


wall 


where, S R = 




<25 


7 T ’*« S25 

k r 


( 20 ) 


and K r is the equivalent sand grain roughness height in turbulent wall units. This 
condition on co is attributable to Wilcox. A discussion on the physical reasoning behind the 
condition can be found in Wilcox (1988). 

An upper limit is imposed on the value of co at the wall using the following boundary 
condition suggested by Menter(1993) and found effective by Chima(1996), 


= |0 6 _v_ 
W Re P' A / 


( 21 ) 

( 22 ) 

(23) 
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2.5. Time Stepping Scheme 

The discretized governing equations are marched in time to steady state using a Runge- 
Kutta time stepping scheme. Typically a four or five stage scheme is used. To accelerate 
convergence to steady state, the higher order time accuracy of the schemes is sacrificed in 
favor of near optimal high-frequency damping properties of the scheme, which improves 
effectiveness of multi-grid schemes without affecting the accuracy of the final steady state 
solution. The schemes used here were originally proposed for the Euler equations by Jameson, 
Schmidt and Turkel (1981) and have been widely used. 

2.6. Convergence Acceleration 

2.6.1. Local Time Stepping 

For steady flows, Ar is viewed simply as a relaxation parameter that we can take advantage 
of to accelerate the progression of the computed solution to steady state. Thus, rather than use 
the same time step for every equation, time accuracy is sacrificed and a “local” time step is 
used. For non-dissipative, hyperbolic systems of equations, steady state is reached once all 
unsteady waves have been expelled out of the domain. Thus, the local time step is chosen 
based on the local stability limit which, for such systems, is proportional to the time it takes a 
characteristic wave to travel across a cell, i.e., the local time step is computed from constant 
multiple of the local CFL number. For the Navier-Stokes equations, the presence of viscous 
terms lead to more complicated stability limits so that empirical formulas are used, (see 
Jameson, etal., 1981, Amone, et al., 1990, and Swanson and Turkel, 1987) 

For a given grid system, higher flow velocities lead to larger eigenvalues and, therefore, to 
smaller allowable time step size. Higher Reynolds number also reduces the allowable time 
step size. In boundary layers with highly clustered grids, the viscous time step limits the 
allowable time step size. 

2.6.2. Implicit Residual Smoothing 

To further improve stability characteristics and rate of convergence to steady state, a 
technique known as implicit residual smoothing (IRS) is applied in every stage of the Runge- 
Kutta time stepping scheme. Implicit Residual Smoothing was first introduced by Lerat 
(1979) for Lax-Wendrof-type of schemes and by Jameson (1983a,b) for multi-stage schemes. 
Modifications to the scheme by Jameson were introduced by Swanson and Turkel (1987) and 
by Martinelli and Jameson (1988) to improve the performance of the schemes on high aspect- 
ratio grids such as used in simulations of viscous flows. Here the scheme of Martinelli and 
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Jameson is used. For an in-depth discussion on implicit residual smoothing, see Swanson and 
Turkel (1997). 

For multi-block grid systems, the Implicit Residual Smoothing is applied independently 
within each block. For relatively few blocks, this is certainly adequate and does not hamper 
performance (see e.g., Steinthorsson, et al. 1993). For a large number of blocks, performance 
gain due to IRS may be affected, in particular if the blocks are small. However, the cost of 
coupling the blocks in the IRS will likely offset any gain in performance. 

Finally, note that the boundary conditions on the smoothed residual are needed. Here the 
systems of equations are solved assuming that the smoothed residual in the ghost cells is zero. 
Several other approaches have been tried, including implicit treatment of boundary 
conditions. These approaches have been found to provide benefit in some cases but to be 
detrimental in others. Overall, the simple boundary treatment currently used has proven to be 
the most robust. 

2.6.3. Multigrid 

Multigrid convergence acceleration involves the use of a sequence of coarser and coarser 
grids to accelerate the evolution of the computed fine-grid solution to steady state. The 
multigrid concept can be traced at least back to the work of Fedorenko (1964). Multigrid was 
introduced to Computational Fluid Mechanics primarily through the work of Brandt (see 
Brandt, 1977), Jameson, Ni, (see Ni, 1981) and Turkel. The multigrid scheme used in the 
current code is based on the work of Jameson (1983) and was formulated by Amone and 
Swanson (1988) in the original TRAF2D code. 

Although the available theory behind multigrid methods will not be reviewed in detail, it is 
useful for readers unfamiliar with the method to review the underlying ideas. For an 
elementary but enlightening discussion, the interested reader is directed to a tutorial by Briggs 
(1987). For a more advanced discussion of multigrid method applied to Computational Fluid 
Dynamics, see e.g., Hacbusch (1982). 

3. Multi-Block Grid Generation Techniques 

For the numerical solution of problems with complex geometry, grid generation becomes 
a very high priority. So important is the role of grid generation that it dictates the type of 
solver to be used. Currently, the two most popular types of solvers are structured multi-block 
and unstructured tetrahedron. Other techniques include overlapping structured grids, known 
as Chimera grids (see Steger); as well as a combination of structured and unstructured grids. 
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i.e. DRAGON grids (see Kao and Liou). If the question of grid generation is set aside, the 
structured multi-block approach is the clear choice in terms of numerical efficiency and 
memory usage. Regardless of how fast or inexpensive computers become, this advantage will 
remain. 

The generation of unstructured tetrahedron grids is currently regarded as requiring very 
few man hours relative to other types of grids. However, for grids which require anisotropic 
grid resolution, the complexity in grid generation remains. The generation of tetrahedrons 
with aspect ratio far from unity is a developing technology and not widely available. In 
addition, the accuracy of solvers on these types of grids is not well documented. In particular, 
the prediction of skin friction and heat transfer may require significantly more resolution than 
a hexahedron mesh. Ideas such as generating prismatic meshes near viscous boundaries may 
improve accuracy, but will further complicate the grid generation process. 

One could assume that structured multi-block grid generation will remain too tedious and 
labor intensive. It is much more likely, however, that fully automatic structured multi-block 
grid generation will be realized in the near future. The most difficult task to automate is the 
topology generation, that is how the blocks should be arranged and connected to best resolve 
the features of the problem. If one imagines that each block is a cell in an unstructured 
hexahedron mesh, then it becomes apparent that unstructured hexahedron mesh generators 
such as CUBIT or HOUDINI, which are primarily used for finite element analysis in solids, 
may be exploited to generate the topology which could then be used in an existing structured 
multi-block grid generation programs such as GRIDGEN or GridPro. 

Currently, the GridPro software from Program Development Corporation is used 
extensively at NASA Glenn for turbine related calculations. Each new geometry requires the 
development of an appropriate topology. Some general guidelines which help lead to high 
quality grids are as follows: 

• Limit topological singularities to one away from optimal. For instance, in 2D four blocks 
is optimal so having three or five would be acceptable. 

• Have grid lines follow the surface when ever possible. This is especially important for vis- 
cous grids, to keep small grid spacings from being introduced into the main flow path 
where they are not needed. 

• Keep topological singularities off boundaries to maximize grid quality where numerical 
accuracy is most needed. 
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4. Block Merging with the Method of Weakest Descent 

4.1. Introduction 

A methodology for automatic block merging is developed to the point where a 
complicated multi-block grid system is supplied and a merged system with connectivity 
information is returned. The objective of the method is to produce the minimum number of 
blocks. The recently developed Method of Weakest Descent (MWD) is described and its 
application to representative test cases is presented, Rigby (1996) and Rigby et al (1997). The 
MWD is based on the premise that the more internal faces of a multi-block grid that can be 
removed, the lower the final number of blocks. With each internal face that is removed, other 
faces are disqualified from being removed as a result of the merging of two blocks. The MWD 
chooses which internal face to remove so as to minimize the number of disqualified faces, thus 
diminishing the set of available faces at the slowest rate (i.e. weakest descent). Reducing the 
available faces at the slowest rate allows more blocks to be merged before all valid internal 
faces are removed. When more than one candidate disqualifies the same number of other 
faces, the choice among those candidates is made randomly. Because of this randomness, each 
application of the MWD can produce different results. Many tests have been done, with the 
number of initial blocks ranging from 12 to 7936. The tests have shown that, for most cases, 
every application of the MWD produces a result near the expected minimum. Repeated 
application increases the likelihood of realizing the actual minimum. 

4.2. Description of the Method of Weakest Descent 

The MWD procedure is fully three dimensional and uses only the connectivity data so it 
requires very little memory. The method progresses by successively merging blocks which 
result in logically rectangular blocks. No invalid blocks are formed at any stage of the merging 
procedure. The patch which joins two blocks that can be merged is referred to as a candidate. 
The removal of a candidate implies the merging of the blocks on either side of it. In general, 
the objective is to deplete the set of candidates available for removal at the slowest rate 
possible, thus allowing for the most blocks to be merged. It is this notion of depleting the set 
of candidates available for removal at the slowest rate possible which led to the name of the 
algorithm. It is important to realize that when two blocks are merged not only is the candidate 
between the blocks lost, but also several other candidates may be disqualified. The key is to, at 
each stage, merge the blocks which disqualify the fewest candidates. The discussion of the 
simple circle in a box case, in the results section, covers this point in more detail. 
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A flow chart for the algorithm is shown in figure 1. The thick line with filled arrows shows 
the main loop which defines the algorithm. Each execution of this loop is referred to as a stage 
and represents the reduction in the number of blocks by one. 

Each stage of the procedure consists of the following seven steps: 

1. Join any group of patches which cover the entire face of both neighboring 
blocks. 

2. Update patch array which relates patches to blocks, and block faces. With the 
connectivity data, given the patch number, one finds the blocks and associated 
indices. The patch array returns the patch number given the block and face 
information. 

3. From the set of internal faces, determine those which are candidates for 
removal. Several straight forward tests are done to ensure that merging two 
blocks will result in a logically rectangular block. One of the less obvious tests 
is to make sure that a patch which has the same block on both sides, as in an O- 
type grid, is not taken to be a candidate. 

4. If there are any candidates, continue. 

5. Grade each candidate for removal based on how many other candidates would 
be disqualified by its removal. The grade given a candidate is simply the differ- 
ence between the current number of candidates and the number that would exist 
if the two blocks were merged. It is important to include the production of can- 
didates, which may result after groups of patches are joined, in the grade. Also, 
the loss of a candidate due to the formation of an O-type grid is not included in 
the grade to encourage the formation of O-type grids. 

6. Remove a low grade face and merge the neighboring blocks. If there is more 
than one candidate with the same low grade, then one of them is chosen ran- 
domly. 

7. Return to step 1 and repeat until no candidates for removal are found. 

In the present work, grading of candidates was based solely on how many other candidates 
would be disqualified. Additional grading criterion could also be added such as: setting an 
upper limit on block size to meet a memory constraint, or setting a lower limit on the number 
of blocks based on number of available CPUs, or consider the size of the resulting block 


14 


Rigby 


relative to other blocks to enhance load balancing. Then, depending on the reliability of the 
grading procedure, the candidate with the lowest score may be chosen or the random choice 
could be made from a range of scores. In fact, the random choice could be designed in such a 
way as to favor the candidates with lower grades. These ideas will be explored in future work. 

4.3. Sample Applications of the MWD 

4.3.1. Overview 

First, a very simple example of a circle in a box is shown to explain the progression of an 
application of the Method of Weakest Descent. Then, several complex examples are shown. 

For each of the three complex examples, the topology is shown to demonstrate the 
complexity. The effectiveness of the method is presented in terms of the probability versus the 
final number of blocks (where the probability refers to the probability of realizing a particular 
final number of blocks with each application of the MWD). This probability is calculated by 
performing 500 realizations of the MWD; then, the number of times each final number of 
blocks is realized is divided by 500. 

To further demonstrate the effectiveness of the MWD, two other grading procedures are 
also presented. First, the MWD was run as previously described. The grading procedure was 
then disabled with the candidate for removal chosen randomly from the entire set of valid 
internal faces at each stage. Finally, instead of choosing the low grade, the high grade 
candidate was chosen. 

A discussion on the CPU and memory requirements is also presented with regard to the 
third complex example. 

4.3.2. Circle in a Box 

Figure 2 shows a topology for a circle in a box. This case is useful to demonstrate the step 
by step progression of the MWD. The solid lines in figure 2 show the outline of blocks that 
will produce a high quality grid for this problem. Near the fiat walls the grid will be nearly 
cartesian, while near the circle the grid will be nearly polar. 

Figure 3 shows the state of the blocks at each stage of the MWD as the grid is merged 
from the initial twelve blocks to the minimum of five blocks. The number in the gray arrow to 
the left of each configuration refers to the stage. The initial configuration is shown in the upper 
left of figure 3. Initially, each patch is a candidate and receives a grade of either 3 or 5. One of 
the candidates with grade 3 is chosen randomly and is removed, as shown at stage 1 in figure 
3. In stage 1, the circled candidate from the initial configuration has been removed, and two 
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other patches which were previously candidates were disqualified as indicated by the x’s. The 
grades of patches in neighboring blocks has also been updated to reflect the change in the 
topology. The progression to stage 2 is accomplished by randomly choosing one of the 
patches with a grade of 2. Removal of candidates continues until stage 7, when no more 
candidates exist. The final number of blocks is the initial number of blocks minus the number 
of stages, that is 5. 

Figure 4 shows the number of candidates available at each stage drawn in figure 3. At each 
stage, the slope of this curve is equal to the negative of the grade given the candidate removed, 
except for stage 6 where a candidate of grade 2 was removed, but a loss of 3 was incurred. The 
reason for this is because, as mentioned earlier, the loss of a candidate due to the formation of 
an O-type grid is not included in the grade. The purpose of figure 4 is to clearly demonstrate 
the need to minimize the rate at which candidates are depleted so as to allow the maximum 
number of stages. Recall that each stage marks the reduction in the number of blocks by one. 

If at each stage the highest grade candidate were chosen for removal, one would find that 
the final number of blocks would be eight, since only four stages would be possible. 

4.3.3. Internal Passage with Ribs and Bleed 

The final example is a complex internal coolant passage with ribs and bleed holes. Figure 
9a shows the surface grid constructed using the GridPro software[l]. Notice that, by taking 
full advantage of the multi-block capability, high quality grids can be formed around the ribs 
and holes. Figure 9b shows the block boundaries for the initial 282 block grid from GridPro. 
Figure 9c shows the block boundaries of the merged grid system, which has only 22 blocks. 

Figure 10 shows the probability results for 500 realizations. The results are very 
encouraging. The likelihood of producing a result near the observed minimum is very high 
with every application of the MWD. 

Figure 1 1 shows the number of candidates available at each stage versus the number of 
blocks at each stage. Ten realizations are shown for each grading procedure. Clearly, the 
MWD produces the slowest loss of candidates, thus allowing more blocks to be merged before 
the pool of candidates is depleted. 

To investigate the CPU time requirements as the number of blocks increases to much 
higher levels, the topology for a case similar to this one was divided into additional blocks. 
The original topology had 124 blocks. The divided topologies had 992, 1776, and 7936 
blocks. The original topology had 238 internal patches, while the divided topologies had 
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2440, 4700, and 21664, respectively. Figure 12 shows the CPU time required for a single 
application of the MWD for each of the divided topologies on three different machines. The 
tests were run on a Silicon Graphics Indigo 2 with an R4000 chip, a Silicon Graphics Power 
Challenge with an R10000 chip, and on the NAS CRAY C90 (eagle). From figure 11, it can be 
seen that on the order of 5 seconds is required to merge 1000 blocks with the R4000 chip or on 
the C90. The R10000 chip consistently performs five to six times faster than the other two 
machines. The CPU time is seen to increase like the square of the number of blocks. It is 
expected that improvements in the implementation of the algorithm should allow the CPU 
time to increase linearly with the number of blocks. 

The reason that the C90 does not show an advantage over the Silicon Graphics machines is 
because the code is dominated by operations which are not vectorized. The case with 7936 
blocks was only run on the C90 because of the memory requirements. In its current state, the 
memory requirement for the code increases like the number of blocks raised to the 1.6 power. 
For 1776 block case, 3.2 MW were required, while for the 7936 block case, 35 MW were 
required. It is expected that improvements in the implementation of the algorithm should 
allow the memory to increase linearly with the number of blocks. 

4.3.4. Final Remarks 

A methodology for automatically merging complicated grid systems has been described. 
The method is currently developed to the point where 9 grid system is input, and a merged 
system with the associated connectivity information is produced; no additional user 
intervention is required. The Method of Weakest Descent has been demonstrated as an 
effective automatic block merging procedure. A single application has been shown to very 
quickly produce a result near the expected optimum. With additional applications, results even 
nearer the optimum can result. It is recommended that a specified number of realizations be 
done for each case, with the specified number perhaps being dependent upon the initial, 
number of internal patches. 

Without the automatic block merging capability, topologies as complicated as those 
presented would probably not even be attempted. Calculations on topologies with hundreds of 
blocks would be tedious to set up and post process, and be numerically inefficient to run. The 
ability to automatically merge grid systems allows the person generating the grids to take full 
advantage of the multi-block methodology. The Method of Weakest Descent has been 
incorporated into the GridPro software and is routinely used for turbine related calculations 
done at the NASA Glenn Research Center. 
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BEGIN 


Set number of realizations to be taken 


Read connectivity and block sizes 


Deduce how many patches are on each 

FACE OF EACH BLOCK 


Create an array which relates a patch 

TO THE NTH PATCH ON A GIVEN FACE OF A 
GIVEN BLOCK 


1 . Join any patches that can be joined 


2. Update patch array 


3. Determine which patches are 
CANDIDATES FOR REMOVAL 


4. ARE THERE ANY CANDIDATES? 

v-tYES 


5. Grade candidates on how many 

CANDIDATES WILL BE LOST IF REMOVED 


6. Remove a low grade candidate by 
JOINING ITS NEIGHBORING BLOCKS 


IS THIS THE BEST SCHEDULE SO FAR? 

„ ,YES 

Output merging schedule and 

CONNECTIVITY 


Should another realization be done? 


STOP 


Figure 1. Flow chart for the Method of Weakest Descent. Thick line shows main loop. 
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Figure 3. Stages for a circle in a box case. Grade shown near candidate. Invalid patches have x. Circle shows candi 
date which is about to be removed. 
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Figure 4. Number of available candidates at each stage 
for the circle in a box example. 
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